<!-- HTML header for doxygen 1.8.13-->
<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN" "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd">
<html xmlns="http://www.w3.org/1999/xhtml">
<head>
<meta http-equiv="Content-Type" content="text/xhtml;charset=UTF-8"/>
<meta http-equiv="X-UA-Compatible" content="IE=9"/>
<meta name="generator" content="Doxygen 1.8.14"/>
<meta name="viewport" content="width=device-width, initial-scale=1"/>
<title>Taskflow Handbook</title>
<link href="tabs.css" rel="stylesheet" type="text/css"/>
<link rel="icon" type="image/x-icon" href="favicon.ico" />
<script type="text/javascript" src="jquery.js"></script>
<script type="text/javascript" src="dynsections.js"></script>
<link href="navtree.css" rel="stylesheet" type="text/css"/>
<script type="text/javascript" src="resize.js"></script>
<script type="text/javascript" src="navtreedata.js"></script>
<script type="text/javascript" src="navtree.js"></script>
<script type="text/javascript">
/* @license magnet:?xt=urn:btih:cf05388f2679ee054f2beb29a391d25f4e673ac3&amp;dn=gpl-2.0.txt GPL-v2 */
  $(document).ready(initResizable);
/* @license-end */</script>
<link href="search/search.css" rel="stylesheet" type="text/css"/>
<script type="text/javascript" src="search/searchdata.js"></script>
<script type="text/javascript" src="search/search.js"></script>
<link href="doxygen.css" rel="stylesheet" type="text/css" />
</head>
<body>
<div id="top"><!-- do not remove this div, it is closed by doxygen! -->
<div id="titlearea">
<table cellspacing="0" cellpadding="0">
 <tbody>
 <tr style="height: 56px;">
  <td id="projectalign" style="padding-left: 0.5em;">
   <div id="projectname"><a href="https://taskflow.github.io/">Taskflow</a>
   &#160;<span id="projectnumber">3.0.0-Master-Branch</span>
   </div>
  </td>
 </tr>
 </tbody>
</table>
</div>
<!-- end header part -->
<!-- Generated by Doxygen 1.8.14 -->
<script type="text/javascript">
/* @license magnet:?xt=urn:btih:cf05388f2679ee054f2beb29a391d25f4e673ac3&amp;dn=gpl-2.0.txt GPL-v2 */
var searchBox = new SearchBox("searchBox", "search",false,'Search');
/* @license-end */
</script>
<script type="text/javascript" src="menudata.js"></script>
<script type="text/javascript" src="menu.js"></script>
<script type="text/javascript">
/* @license magnet:?xt=urn:btih:cf05388f2679ee054f2beb29a391d25f4e673ac3&amp;dn=gpl-2.0.txt GPL-v2 */
$(function() {
  initMenu('',true,false,'search.php','Search');
  $(document).ready(function() { init_search(); });
});
/* @license-end */</script>
<div id="main-nav"></div>
</div><!-- top -->
<div id="side-nav" class="ui-resizable side-nav-resizable">
  <div id="nav-tree">
    <div id="nav-tree-contents">
      <div id="nav-sync" class="sync"></div>
    </div>
  </div>
  <div id="splitbar" style="-moz-user-select:none;" 
       class="ui-resizable-handle">
  </div>
</div>
<script type="text/javascript">
/* @license magnet:?xt=urn:btih:cf05388f2679ee054f2beb29a391d25f4e673ac3&amp;dn=gpl-2.0.txt GPL-v2 */
$(document).ready(function(){initNavTree('kmeans.html','');});
/* @license-end */
</script>
<div id="doc-content">
<!-- window showing the filter options -->
<div id="MSearchSelectWindow"
     onmouseover="return searchBox.OnSearchSelectShow()"
     onmouseout="return searchBox.OnSearchSelectHide()"
     onkeydown="return searchBox.OnSearchSelectKey(event)">
</div>

<!-- iframe showing the search results (closed by default) -->
<div id="MSearchResultsWindow">
<iframe src="javascript:void(0)" frameborder="0" 
        name="MSearchResults" id="MSearchResults">
</iframe>
</div>

<div class="header">
  <div class="headertitle">
<div class="title">k-means Clustering </div>  </div>
</div><!--header-->
<div class="contents">
<div class="textblock"><p>We study a fundamental clustering problem in unsupervised learning, <em>k-means clustering</em>. We will begin by discussing the problem formulation and then learn how to write a parallel k-means algorithm.</p>
<h1><a class="anchor" id="KMeansProblemFormulation"></a>
Problem Formulation</h1>
<p>k-means clustering uses <em>centroids</em>, k different randomly-initiated points in the data, and assigns every data point to the nearest centroid. After every point has been assigned, the centroid is moved to the average of all of the points assigned to it. We describe the k-means algorithm in the following steps:</p>
<ul>
<li>
Step 1: initialize k random centroids </li>
<li>
Step 2: for every data point, find the nearest centroid (L2 distance or other measurements) and assign the point to it </li>
<li>
Step 3: for every centroid, move the centroid to the average of the points assigned to that centroid </li>
<li>
Step 4: go to Step 2 until converged (no more changes in the last few iterations) or maximum iterations reached </li>
</ul>
<p>The algorithm is illustrated as follows:</p>
<div class="image">
<img src="kmeans_1.png" alt="kmeans_1.png" width="100%"/>
</div>
<p>A sequential implementation of k-means is described as follows:</p>
<div class="fragment"><div class="line"><span class="comment">// sequential implementation of k-means on a CPU</span></div><div class="line"><span class="comment">// N: number of points, K: number of clusters, M: number of iterations, px/py: 2D point vector </span></div><div class="line"><span class="keywordtype">void</span> kmeans_seq(<span class="keywordtype">int</span> N, <span class="keywordtype">int</span> K, <span class="keywordtype">int</span> M, <span class="keyword">const</span> <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;float&gt;</a>&amp; px,  <span class="keyword">const</span> <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;float&gt;</a>&amp; py) {</div><div class="line"></div><div class="line">  <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;int&gt;</a> c(K);</div><div class="line">  <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;float&gt;</a> sx(K), sy(K), mx(K), my(K);</div><div class="line"></div><div class="line">  <span class="comment">// initial centroids</span></div><div class="line">  <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/copy_n.html">std::copy_n</a>(px.begin(), K, mx.begin());</div><div class="line">  <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/copy_n.html">std::copy_n</a>(py.begin(), K, my.begin());</div><div class="line">  </div><div class="line">  <span class="comment">// k-means iteration</span></div><div class="line">  <span class="keywordflow">for</span>(<span class="keywordtype">int</span> m=0; m&lt;M; m++) {</div><div class="line"></div><div class="line">    <span class="comment">// clear the storage</span></div><div class="line">    <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/fill_n.html">std::fill_n</a>(sx.begin(), K, 0.0f);</div><div class="line">    <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/fill_n.html">std::fill_n</a>(sy.begin(), K, 0.0f);</div><div class="line">    <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/fill_n.html">std::fill_n</a>(c.begin(), K, 0);</div><div class="line"></div><div class="line">    <span class="comment">// find the best k (cluster id) for each point</span></div><div class="line">    <span class="keywordflow">for</span>(<span class="keywordtype">int</span> i=0; i&lt;N; ++i) {</div><div class="line">      <span class="keywordtype">float</span> x = px[i];</div><div class="line">      <span class="keywordtype">float</span> y = py[i];</div><div class="line">      <span class="keywordtype">float</span> best_d = <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/types/numeric_limits/max.html">std::numeric_limits&lt;float&gt;::max</a>();</div><div class="line">      <span class="keywordtype">int</span> best_k = 0;</div><div class="line">      <span class="keywordflow">for</span> (<span class="keywordtype">int</span> k = 0; k &lt; K; ++k) {</div><div class="line">        <span class="keyword">const</span> <span class="keywordtype">float</span> d = L2(x, y, mx[k], my[k]);</div><div class="line">        <span class="keywordflow">if</span> (d &lt; best_d) {</div><div class="line">          best_d = d;</div><div class="line">          best_k = k;</div><div class="line">        }</div><div class="line">      }</div><div class="line">      sx[best_k] += x;</div><div class="line">      sy[best_k] += y;</div><div class="line">      c [best_k] += 1;</div><div class="line">    }</div><div class="line"></div><div class="line">    <span class="comment">// update the centroid</span></div><div class="line">    <span class="keywordflow">for</span>(<span class="keywordtype">int</span> k=0; k&lt;K; k++) {</div><div class="line">      <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a> = <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/max.html">max</a>(1, c[k]);  <span class="comment">// turn 0/0 to 0/1</span></div><div class="line">      mx[k] = sx[k] / <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a>;</div><div class="line">      my[k] = sy[k] / <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a>;</div><div class="line">    }</div><div class="line">  }</div><div class="line"></div><div class="line">  <span class="comment">// print the k centroids found</span></div><div class="line">  <span class="keywordflow">for</span>(<span class="keywordtype">int</span> k=0; k&lt;K; ++k) {</div><div class="line">    <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/io/basic_ostream.html">std::cout</a> &lt;&lt; <span class="stringliteral">&quot;centroid &quot;</span> &lt;&lt; k &lt;&lt; <span class="stringliteral">&quot;: &quot;</span> &lt;&lt; <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/io/manip/setw.html">std::setw</a>(10) &lt;&lt; mx[k] &lt;&lt; <span class="charliteral">&#39; &#39;</span></div><div class="line">                                          &lt;&lt; <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/io/manip/setw.html">std::setw</a>(10) &lt;&lt; my[k] &lt;&lt; <span class="charliteral">&#39;\n&#39;</span>;</div><div class="line">  }</div><div class="line">}</div></div><!-- fragment --><h1><a class="anchor" id="ParallelKMeansUsingCPUs"></a>
Parallel k-means using CPUs</h1>
<p>The second step of k-means algorithm, <em>assigning every point to the nearest centroid</em>, is highly parallelizable across individual points. We can create a <em>parallel-for</em> task to run parallel iterations.</p>
<div class="fragment"><div class="line"><a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;int&gt;</a> best_ks(N);  <span class="comment">// nearest centroid of each point</span></div><div class="line"></div><div class="line"><span class="keywordtype">unsigned</span> P = 12;  <span class="comment">// 12 partitioned tasks</span></div><div class="line"></div><div class="line"><span class="comment">// update cluster</span></div><div class="line">taskflow.parallel_for(0, N, 1, [&amp;](<span class="keywordtype">int</span> i){</div><div class="line">  <span class="keywordtype">float</span> x = px[i];</div><div class="line">  <span class="keywordtype">float</span> y = py[i];</div><div class="line">  <span class="keywordtype">float</span> best_d = <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/types/numeric_limits/max.html">std::numeric_limits&lt;float&gt;::max</a>();</div><div class="line">  <span class="keywordtype">int</span> best_k = 0;</div><div class="line">  <span class="keywordflow">for</span> (<span class="keywordtype">int</span> k = 0; k &lt; K; ++k) {</div><div class="line">    <span class="keyword">const</span> <span class="keywordtype">float</span> d = L2(x, y, mx[k], my[k]);</div><div class="line">    <span class="keywordflow">if</span> (d &lt; best_d) {</div><div class="line">      best_d = d;</div><div class="line">      best_k = k;</div><div class="line">    }</div><div class="line">  }</div><div class="line">  best_ks[i] = best_k;</div><div class="line">});</div></div><!-- fragment --><p>The third step of moving every centroid to the average of points is also parallelizable across individual centroids. However, since k is typically not large, one task of doing this update is sufficient.</p>
<div class="fragment"><div class="line">taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){</div><div class="line">  <span class="comment">// sum of points</span></div><div class="line">  <span class="keywordflow">for</span>(<span class="keywordtype">int</span> i=0; i&lt;N; i++) {</div><div class="line">    sx[best_ks[i]] += px[i];</div><div class="line">    sy[best_ks[i]] += py[i];</div><div class="line">    c [best_ks[i]] += 1;</div><div class="line">  }</div><div class="line">  </div><div class="line">  <span class="comment">// average of points</span></div><div class="line">  <span class="keywordflow">for</span>(<span class="keywordtype">int</span> k=0; k&lt;K; ++k) {</div><div class="line">    <span class="keyword">auto</span> <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a> = <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/max.html">max</a>(1, c[k]);  <span class="comment">// turn 0/0 to 0/1</span></div><div class="line">    mx[k] = sx[k] / <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a>;</div><div class="line">    my[k] = sy[k] / <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a>;</div><div class="line">  }</div><div class="line">});</div></div><!-- fragment --><p>To describe <code>M</code> iterations, we create a condition task that loops the second step of the algorithm by <code>M</code> times. The return value of zero goes to the first successor which we will connect to the task of the second step later; otherwise, k-means completes.</p>
<div class="fragment"><div class="line">taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([m=0, M]() <span class="keyword">mutable</span> {</div><div class="line">  <span class="keywordflow">return</span> (m++ &lt; M) ? 0 : 1;</div><div class="line">});</div></div><!-- fragment --><p>The entire code of CPU-parallel k-means is shown below. Here we use an additional storage, <code>best_ks</code>, to record the nearest centroid of a point at an iteration.</p>
<div class="fragment"><div class="line"><span class="comment">// N: number of points, K: number of clusters, M: number of iterations, px/py: 2D point vector </span></div><div class="line"><span class="keywordtype">void</span> kmeans_par(<span class="keywordtype">int</span> N, <span class="keywordtype">int</span> K, <span class="keywordtype">int</span> M, cconst <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;float&gt;</a>&amp; px,  <span class="keyword">const</span> <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;float&gt;</a>&amp; py) {</div><div class="line"></div><div class="line">  <span class="keywordtype">unsigned</span> P = 12;  <span class="comment">// 12 partitions of the parallel-for graph</span></div><div class="line"></div><div class="line">  <a class="code" href="classtf_1_1Executor.html">tf::Executor</a> executor;</div><div class="line">  <a class="code" href="classtf_1_1Taskflow.html">tf::Taskflow</a> taskflow(<span class="stringliteral">&quot;K-Means&quot;</span>);</div><div class="line"></div><div class="line">  <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;int&gt;</a> c(K), best_ks(N);</div><div class="line">  <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;float&gt;</a> sx(K), sy(K), mx(K), my(K);</div><div class="line"></div><div class="line">  <span class="comment">// initial centroids</span></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> init = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){</div><div class="line">    <span class="keywordflow">for</span>(<span class="keywordtype">int</span> i=0; i&lt;K; ++i) {</div><div class="line">      mx[i] = px[i];</div><div class="line">      my[i] = py[i];</div><div class="line">    }</div><div class="line">  }).name(<span class="stringliteral">&quot;init&quot;</span>);</div><div class="line"></div><div class="line">  <span class="comment">// clear the storage</span></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> clean_up = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){</div><div class="line">    <span class="keywordflow">for</span>(<span class="keywordtype">int</span> k=0; k&lt;K; ++k) {</div><div class="line">      sx[k] = 0.0f;</div><div class="line">      sy[k] = 0.0f;</div><div class="line">      c [k] = 0;</div><div class="line">    }</div><div class="line">  }).name(<span class="stringliteral">&quot;clean_up&quot;</span>);</div><div class="line"></div><div class="line">  <span class="comment">// update cluster</span></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> pf = taskflow.parallel_for(0, N, 1, [&amp;](<span class="keywordtype">int</span> i){</div><div class="line">    <span class="keywordtype">float</span> x = px[i];</div><div class="line">    <span class="keywordtype">float</span> y = py[i];</div><div class="line">    <span class="keywordtype">float</span> best_d = <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/types/numeric_limits/max.html">std::numeric_limits&lt;float&gt;::max</a>();</div><div class="line">    <span class="keywordtype">int</span> best_k = 0;</div><div class="line">    <span class="keywordflow">for</span> (<span class="keywordtype">int</span> k = 0; k &lt; K; ++k) {</div><div class="line">      <span class="keyword">const</span> <span class="keywordtype">float</span> d = L2(x, y, mx[k], my[k]);</div><div class="line">      <span class="keywordflow">if</span> (d &lt; best_d) {</div><div class="line">        best_d = d;</div><div class="line">        best_k = k;</div><div class="line">      }</div><div class="line">    }</div><div class="line">    best_ks[i] = best_k;</div><div class="line">  }).name(<span class="stringliteral">&quot;parallel-for&quot;</span>);</div><div class="line"></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> update_cluster = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){</div><div class="line">    <span class="keywordflow">for</span>(<span class="keywordtype">int</span> i=0; i&lt;N; i++) {</div><div class="line">      sx[best_ks[i]] += px[i];</div><div class="line">      sy[best_ks[i]] += py[i];</div><div class="line">      c [best_ks[i]] += 1;</div><div class="line">    }</div><div class="line"></div><div class="line">    <span class="keywordflow">for</span>(<span class="keywordtype">int</span> k=0; k&lt;K; ++k) {</div><div class="line">      <span class="keyword">auto</span> <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a> = <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/max.html">max</a>(1, c[k]);  <span class="comment">// turn 0/0 to 0/1</span></div><div class="line">      mx[k] = sx[k] / <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a>;</div><div class="line">      my[k] = sy[k] / <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a>;</div><div class="line">    }</div><div class="line">  }).name(<span class="stringliteral">&quot;update_cluster&quot;</span>);</div><div class="line">  </div><div class="line">  <span class="comment">// convergence check</span></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> condition = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([m=0, M]() <span class="keyword">mutable</span> {</div><div class="line">    <span class="keywordflow">return</span> (m++ &lt; M) ? 0 : 1;</div><div class="line">  }).name(<span class="stringliteral">&quot;converged?&quot;</span>);</div><div class="line"></div><div class="line">  init.<a class="code" href="classtf_1_1Task.html#a8c78c453295a553c1c016e4062da8588">precede</a>(clean_up);</div><div class="line"></div><div class="line">  clean_up.<a class="code" href="classtf_1_1Task.html#a8c78c453295a553c1c016e4062da8588">precede</a>(pf);</div><div class="line">  pf.<a class="code" href="classtf_1_1Task.html#a8c78c453295a553c1c016e4062da8588">precede</a>(update_cluster);</div><div class="line"></div><div class="line">  condition.<a class="code" href="classtf_1_1Task.html#a8c78c453295a553c1c016e4062da8588">precede</a>(clean_up)</div><div class="line">           .<a class="code" href="classtf_1_1Task.html#a331b1b726555072e7c7d10941257f664">succeed</a>(update_cluster);</div><div class="line"></div><div class="line">  executor.<a class="code" href="classtf_1_1Executor.html#a81f35d5b0a20ac0646447eb80d97c0aa">run</a>(taskflow).wait();</div><div class="line">}</div></div><!-- fragment --><p>The taskflow consists of two parts, a <code>clean_up</code> task and a parallel-for graph. The former cleans up the storage <code>sx</code>, <code>sy</code>, and <code>c</code> that are used to average points for new centroids, and the later parallelizes the searching for nearest centroids across individual points using 12 tasks (may vary depending on the machine). If the iteration count is smaller than <code>M</code>, the condition task returns 0 to let the execution path go back to <code>clean_up</code>. Otherwise, it returns 1 to stop (i.e., no successor tasks at index 1). The taskflow graph is illustrated below:</p>
<div class="image">
<object type="image/svg+xml" data="kmeans_2.svg" width="100%">kmeans_2.svg</object>
</div>
<p>The scheduler starts with <code>init</code>, moves on to <code>clean_up</code>, and then enters the parallel-for task <code>paralle-for</code> that spawns a subflow of 12 workers to perform parallel iterations. When <code>parallel-for</code> completes, it updates the cluster centroids and checks if they have converged through a condition task. If not, the condition task informs the scheduler to go back to <code>clean_up</code> and then <code>parallel-for</code>; otherwise, it returns a nominal index to stop the scheduler.</p>
<h1><a class="anchor" id="ParallelKMeansUsingGPUs"></a>
Parallel k-means using GPUs</h1>
<p>We observe Step 2 and Step 3 of the algorithm are parallelizable across individual points for use to harness the power of GPU:</p>
<ol>
<li>
for every data point, find the nearest centroid (L2 distance or other measurements) and assign the point to it </li>
<li>
for every centroid, move the centroid to the average of the points assigned to that centroid. </li>
</ol>
<p>At a fine-grained level, we request one GPU thread to work on one point for Step 2 and one GPU thread to work on one centroid for Step 3.</p>
<div class="fragment"><div class="line"><span class="comment">// px/py: 2D points, N: number of points, mx/my: centroids, K: number of clusters</span></div><div class="line"><span class="comment">// sx/sy/c: storage to compute the average</span></div><div class="line">__global__ <span class="keywordtype">void</span> assign_clusters(</div><div class="line">  <span class="keywordtype">float</span>* px, <span class="keywordtype">float</span>* py, <span class="keywordtype">int</span> N, <span class="keywordtype">float</span>* mx, <span class="keywordtype">float</span>* my, <span class="keywordtype">float</span>* sx, <span class="keywordtype">float</span>* sy, <span class="keywordtype">int</span> K, <span class="keywordtype">int</span>* c</div><div class="line">) {</div><div class="line">  <span class="keyword">const</span> <span class="keywordtype">int</span> index = blockIdx.x * blockDim.x + threadIdx.x;</div><div class="line"></div><div class="line">  <span class="keywordflow">if</span> (index &gt;= N) {</div><div class="line">    <span class="keywordflow">return</span>;</div><div class="line">  }</div><div class="line"></div><div class="line">  <span class="comment">// Make global loads once.</span></div><div class="line">  <span class="keywordtype">float</span> x = px[index];</div><div class="line">  <span class="keywordtype">float</span> y = py[index];</div><div class="line"></div><div class="line">  <span class="keywordtype">float</span> best_dance = FLT_MAX;</div><div class="line">  <span class="keywordtype">int</span> best_k = 0;</div><div class="line">  <span class="keywordflow">for</span> (<span class="keywordtype">int</span> k = 0; k &lt; K; ++k) {</div><div class="line">    <span class="keywordtype">float</span> d = L2(x, y, mx[k], my[k]);</div><div class="line">    <span class="keywordflow">if</span> (d &lt; best_d) {</div><div class="line">      best_d = d;</div><div class="line">      best_k = k;</div><div class="line">    }   </div><div class="line">  }</div><div class="line"></div><div class="line">  atomicAdd(&amp;sx[best_k], x); </div><div class="line">  atomicAdd(&amp;sy[best_k], y); </div><div class="line">  atomicAdd(&amp;c [best_k], 1); </div><div class="line">}</div><div class="line"></div><div class="line"><span class="comment">// mx/my: centroids, sx/sy/c: storage to compute the average</span></div><div class="line">__global__ <span class="keywordtype">void</span> compute_new_means(<span class="keywordtype">float</span>* mx, <span class="keywordtype">float</span>* my, <span class="keywordtype">float</span>* sx, <span class="keywordtype">float</span>* sy, <span class="keywordtype">int</span>* c</div><div class="line">) {</div><div class="line">  <span class="keywordtype">int</span> k = threadIdx.x;</div><div class="line">  <span class="keywordtype">int</span> <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a> = <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/max.html">max</a>(1, c[k]);  <span class="comment">// turn 0/0 to 0/1</span></div><div class="line">  mx[k] = sx[k] / <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a>;</div><div class="line">  my[k] = sy[k] / <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/algorithm/count.html">count</a>;</div><div class="line">}</div></div><!-- fragment --><p>When we recompute the cluster centroids to be the mean of all points assigned to a particular centroid, multiple GPU threads may access the sum arrays, <code>sx</code> and <code>sy</code>, and the count array, <code>c</code>. To avoid data race, we use a simple <code>atomicAdd</code> method. Based on the two kernels, the entire code of CPU-GPU collaborative tasking is described as follows:</p>
<div class="fragment"><div class="line"><span class="comment">// N: number of points, K: number of clusters, M: number of iterations, px/py: 2D point vector </span></div><div class="line"><span class="keywordtype">void</span> kmeans_gpu(<span class="keywordtype">int</span> N, <span class="keywordtype">int</span> K, <span class="keywordtype">int</span> M, cconst <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;float&gt;</a>&amp; px,  <span class="keyword">const</span> <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;float&gt;</a>&amp; py) {</div><div class="line"></div><div class="line">  <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector&lt;float&gt;</a> h_mx, h_my;</div><div class="line">  <span class="keywordtype">float</span> *d_px, *d_py, *d_mx, *d_my, *d_sx, *d_sy, *d_c;</div><div class="line"></div><div class="line">  <span class="keywordflow">for</span>(<span class="keywordtype">int</span> i=0; i&lt;K; ++i) {</div><div class="line">    h_mx.push_back(h_px[i]);</div><div class="line">    h_my.push_back(h_py[i]);</div><div class="line">  }</div><div class="line"></div><div class="line">  <span class="comment">// create a taskflow graph</span></div><div class="line">  <a class="code" href="classtf_1_1Executor.html">tf::Executor</a> executor;</div><div class="line">  <a class="code" href="classtf_1_1Taskflow.html">tf::Taskflow</a> taskflow(<span class="stringliteral">&quot;K-Means&quot;</span>);</div><div class="line">  </div><div class="line">  <span class="comment">// allocate GPU memory</span></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> allocate_px = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){ cudaMalloc(&amp;d_px, N*<span class="keyword">sizeof</span>(<span class="keywordtype">float</span>)); })</div><div class="line">                             .name(<span class="stringliteral">&quot;allocate_px&quot;</span>);</div><div class="line"></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> allocate_py = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){ cudaMalloc(&amp;d_py, N*<span class="keyword">sizeof</span>(<span class="keywordtype">float</span>)); })</div><div class="line">                             .name(<span class="stringliteral">&quot;allocate_py&quot;</span>);</div><div class="line"></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> allocate_mx = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){ cudaMalloc(&amp;d_mx, K*<span class="keyword">sizeof</span>(<span class="keywordtype">float</span>)); })</div><div class="line">                             .name(<span class="stringliteral">&quot;allocate_mx&quot;</span>);</div><div class="line"></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> allocate_my = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){ cudaMalloc(&amp;d_my, K*<span class="keyword">sizeof</span>(<span class="keywordtype">float</span>)); })</div><div class="line">                             .name(<span class="stringliteral">&quot;allocate_my&quot;</span>);</div><div class="line"></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> allocate_sy = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){ cudaMalloc(&amp;d_sy, K*<span class="keyword">sizeof</span>(<span class="keywordtype">float</span>)); })</div><div class="line">                             .name(<span class="stringliteral">&quot;allocate_sy&quot;</span>);</div><div class="line"></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> allocate_c = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){ cudaMalloc(&amp;d_c, K*<span class="keyword">sizeof</span>(<span class="keywordtype">float</span>)); })</div><div class="line">                            .name(<span class="stringliteral">&quot;allocate_c&quot;</span>);</div><div class="line">  </div><div class="line">  <span class="comment">// transfer data from the host to the GPU</span></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> h2d = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](<a class="code" href="classtf_1_1cudaFlow.html">tf::cudaFlow</a>&amp; cf){</div><div class="line">    cf.<a class="code" href="classtf_1_1cudaFlow.html#af03e04771b655f9e629eb4c22e19b19f">copy</a>(d_px, h_px.data(), N).name(<span class="stringliteral">&quot;h2d_px&quot;</span>);</div><div class="line">    cf.<a class="code" href="classtf_1_1cudaFlow.html#af03e04771b655f9e629eb4c22e19b19f">copy</a>(d_py, h_py.data(), N).name(<span class="stringliteral">&quot;h2d_py&quot;</span>);</div><div class="line">    cf.<a class="code" href="classtf_1_1cudaFlow.html#af03e04771b655f9e629eb4c22e19b19f">copy</a>(d_mx, h_mx.data(), K).name(<span class="stringliteral">&quot;h2d_mx&quot;</span>);</div><div class="line">    cf.<a class="code" href="classtf_1_1cudaFlow.html#af03e04771b655f9e629eb4c22e19b19f">copy</a>(d_my, h_my.data(), K).name(<span class="stringliteral">&quot;h2d_my&quot;</span>);</div><div class="line">  }).name(<span class="stringliteral">&quot;h2d&quot;</span>);</div><div class="line">  </div><div class="line">  <span class="comment">// GPU task graph of the main k-means body</span></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> kmeans = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](<a class="code" href="classtf_1_1cudaFlow.html">tf::cudaFlow</a>&amp; cf){</div><div class="line"></div><div class="line">    <a class="code" href="classtf_1_1cudaTask.html">tf::cudaTask</a> zero_c = cf.<a class="code" href="classtf_1_1cudaFlow.html#a40172fac4464f6d805f75921ea3c2a3b">zero</a>(d_c, K).<a class="code" href="classtf_1_1cudaTask.html#ab81b4f71a44af8d61758524f0c274962">name</a>(<span class="stringliteral">&quot;zero_c&quot;</span>);</div><div class="line">    <a class="code" href="classtf_1_1cudaTask.html">tf::cudaTask</a> zero_sx = cf.<a class="code" href="classtf_1_1cudaFlow.html#a40172fac4464f6d805f75921ea3c2a3b">zero</a>(d_sx, K).<a class="code" href="classtf_1_1cudaTask.html#ab81b4f71a44af8d61758524f0c274962">name</a>(<span class="stringliteral">&quot;zero_sx&quot;</span>);</div><div class="line">    <a class="code" href="classtf_1_1cudaTask.html">tf::cudaTask</a> zero_sy = cf.<a class="code" href="classtf_1_1cudaFlow.html#a40172fac4464f6d805f75921ea3c2a3b">zero</a>(d_sy, K).<a class="code" href="classtf_1_1cudaTask.html#ab81b4f71a44af8d61758524f0c274962">name</a>(<span class="stringliteral">&quot;zero_sy&quot;</span>);</div><div class="line"></div><div class="line">    <a class="code" href="classtf_1_1cudaTask.html">tf::cudaTask</a> cluster = cf.<a class="code" href="classtf_1_1cudaFlow.html#adb731be71bdd436dfb5e36e6213a9a17">kernel</a>(</div><div class="line">      (N+1024-1) / 1024, 1024, 0, assign_clusters, d_px, d_py, N, d_mx, d_my, d_sx, d_sy, K, d_c</div><div class="line">    ).<a class="code" href="classtf_1_1cudaTask.html#ab81b4f71a44af8d61758524f0c274962">name</a>(<span class="stringliteral">&quot;cluster&quot;</span>);</div><div class="line"></div><div class="line">    <a class="code" href="classtf_1_1cudaTask.html">tf::cudaTask</a> new_centroid = cf.<a class="code" href="classtf_1_1cudaFlow.html#adb731be71bdd436dfb5e36e6213a9a17">kernel</a>(</div><div class="line">      1, K, 0, compute_new_means, d_mx, d_my, d_sx, d_sy, d_c</div><div class="line">    ).<a class="code" href="classtf_1_1cudaTask.html#ab81b4f71a44af8d61758524f0c274962">name</a>(<span class="stringliteral">&quot;new_centroid&quot;</span>);</div><div class="line"></div><div class="line">    cluster.<a class="code" href="classtf_1_1cudaTask.html#abdd68287ec4dff4216af34d1db44d1b4">precede</a>(new_centroid)</div><div class="line">           .<a class="code" href="classtf_1_1cudaTask.html#a4a9ca1a34bac47e4c9b04eb4fb2f7775">succeed</a>(zero_c, zero_sx, zero_sy);</div><div class="line">  }).name(<span class="stringliteral">&quot;update_means&quot;</span>);</div><div class="line">  </div><div class="line">  <span class="comment">// condition task to check convergence</span></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> condition = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([i=0, M] () <span class="keyword">mutable</span> {</div><div class="line">    <span class="keywordflow">return</span> i++ &lt; M ? 0 : 1;</div><div class="line">  }).name(<span class="stringliteral">&quot;converged?&quot;</span>);</div><div class="line">  </div><div class="line">  <span class="comment">// transfer the result of clusters from GPU to host</span></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> stop = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](<a class="code" href="classtf_1_1cudaFlow.html">tf::cudaFlow</a>&amp; cf){</div><div class="line">    cf.<a class="code" href="classtf_1_1cudaFlow.html#af03e04771b655f9e629eb4c22e19b19f">copy</a>(h_mx.data(), d_mx, K).name(<span class="stringliteral">&quot;d2h_mx&quot;</span>);</div><div class="line">    cf.<a class="code" href="classtf_1_1cudaFlow.html#af03e04771b655f9e629eb4c22e19b19f">copy</a>(h_my.data(), d_my, K).name(<span class="stringliteral">&quot;d2h_my&quot;</span>);</div><div class="line">  }).name(<span class="stringliteral">&quot;d2h&quot;</span>);</div><div class="line">  </div><div class="line">  <span class="comment">// deallocated GPU memory</span></div><div class="line">  <a class="code" href="classtf_1_1Task.html">tf::Task</a> <a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/memory/c/free.html">free</a> = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](){</div><div class="line">    cudaFree(d_px);</div><div class="line">    cudaFree(d_py);</div><div class="line">    cudaFree(d_mx);</div><div class="line">    cudaFree(d_my);</div><div class="line">    cudaFree(d_sx);</div><div class="line">    cudaFree(d_sy);</div><div class="line">    cudaFree(d_c);</div><div class="line">  }).name(<span class="stringliteral">&quot;free&quot;</span>);</div><div class="line"></div><div class="line">  <span class="comment">// build up the dependency</span></div><div class="line">  h2d.<a class="code" href="classtf_1_1Task.html#a331b1b726555072e7c7d10941257f664">succeed</a>(allocate_px, allocate_py, allocate_mx, allocate_my);</div><div class="line"></div><div class="line">  kmeans.<a class="code" href="classtf_1_1Task.html#a331b1b726555072e7c7d10941257f664">succeed</a>(allocate_sx, allocate_sy, allocate_c, h2d)</div><div class="line">        .<a class="code" href="classtf_1_1Task.html#a8c78c453295a553c1c016e4062da8588">precede</a>(condition);</div><div class="line"></div><div class="line">  condition.<a class="code" href="classtf_1_1Task.html#a8c78c453295a553c1c016e4062da8588">precede</a>(kmeans, stop);</div><div class="line"></div><div class="line">  stop.<a class="code" href="classtf_1_1Task.html#a8c78c453295a553c1c016e4062da8588">precede</a>(free);</div><div class="line">  </div><div class="line">  <span class="comment">// dump the taskflow without expanding GPU task graphs</span></div><div class="line">  taskflow.<a class="code" href="classtf_1_1Taskflow.html#a4725d8ea5ff7595d9d71593360538e00">dump</a>(<a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/io/basic_ostream.html">std::cout</a>);</div><div class="line"></div><div class="line">  <span class="comment">// run the taskflow</span></div><div class="line">  executor.<a class="code" href="classtf_1_1Executor.html#a81f35d5b0a20ac0646447eb80d97c0aa">run</a>(taskflow).wait();</div><div class="line">  </div><div class="line">  <span class="comment">// dump the entire taskflow</span></div><div class="line">  taskflow.<a class="code" href="classtf_1_1Taskflow.html#a4725d8ea5ff7595d9d71593360538e00">dump</a>(<a class="codeRef" doxygen="/Users/twhuang/PhD/Code/taskflow/doxygen/cppreference-doxygen-web.tag.xml:http://en.cppreference.com/w/" href="http://en.cppreference.com/w/cpp/io/basic_ostream.html">std::cout</a>);</div><div class="line">}</div></div><!-- fragment --><p>The first dump before executing the taskflow produces the following diagram. The condition tasks introduces a cycle between itself and <code>update_means</code>. Each time it goes back to <code>update_means</code>, the cudaFlow is reconstructed with captured parameters in the closure and offloaded to the GPU.</p>
<div class="image">
<object type="image/svg+xml" data="kmeans_3.svg" width="80%">kmeans_3.svg</object>
</div>
<p>The second dump after executing the taskflow produces the following diagram, with all cudaFlows expanded:</p>
<div class="image">
<object type="image/svg+xml" data="kmeans_4.svg" width="100%">kmeans_4.svg</object>
</div>
<p>The main cudaFlow task, <code>update_means</code>, must not run before all required data has settled down. It precedes a condition task that circles back to itself until we reach <code>M</code> iterations. When iteration completes, the condition task directs the execution path to the cudaFlow, <code>h2d</code>, to copy the results of clusters to <code>h_mx</code> and <code>h_my</code> and then deallocate all GPU memory.</p>
<h1><a class="anchor" id="BuiltInPredicate"></a>
Built-in Predicate</h1>
<p>We observe the GPU task graph parameters remain <em>unchanged</em> across all k-means iterations. In this case, we can leverage <a class="el" href="classtf_1_1cudaFlow.html#aea77b710bf74fb3ccc6043592d4cdbc7" title="offloads the cudaFlow with the given stop predicate and then joins the execution ">tf::cudaFlow::join_until</a> or <a class="el" href="classtf_1_1cudaFlow.html#a9b28ad99e4d3c0208422a2db094df277" title="offloads the cudaFlow and executes it by the given times, and then joins the execution ...">tf::cudaFlow::join_n</a> to run it repeatedly without conditional tasking.</p>
<div class="fragment"><div class="line"><a class="code" href="classtf_1_1Task.html">tf::Task</a> kmeans = taskflow.<a class="code" href="classtf_1_1FlowBuilder.html#a60d7a666cab71ecfa3010b2efb0d6b57">emplace</a>([&amp;](<a class="code" href="classtf_1_1cudaFlow.html">tf::cudaFlow</a>&amp; cf){</div><div class="line"></div><div class="line">  <a class="code" href="classtf_1_1cudaTask.html">tf::cudaTask</a> zero_c = cf.<a class="code" href="classtf_1_1cudaFlow.html#a40172fac4464f6d805f75921ea3c2a3b">zero</a>(d_c, K).<a class="code" href="classtf_1_1cudaTask.html#ab81b4f71a44af8d61758524f0c274962">name</a>(<span class="stringliteral">&quot;zero_c&quot;</span>);</div><div class="line">  <a class="code" href="classtf_1_1cudaTask.html">tf::cudaTask</a> zero_sx = cf.<a class="code" href="classtf_1_1cudaFlow.html#a40172fac4464f6d805f75921ea3c2a3b">zero</a>(d_sx, K).<a class="code" href="classtf_1_1cudaTask.html#ab81b4f71a44af8d61758524f0c274962">name</a>(<span class="stringliteral">&quot;zero_sx&quot;</span>);</div><div class="line">  <a class="code" href="classtf_1_1cudaTask.html">tf::cudaTask</a> zero_sy = cf.<a class="code" href="classtf_1_1cudaFlow.html#a40172fac4464f6d805f75921ea3c2a3b">zero</a>(d_sy, K).<a class="code" href="classtf_1_1cudaTask.html#ab81b4f71a44af8d61758524f0c274962">name</a>(<span class="stringliteral">&quot;zero_sy&quot;</span>);</div><div class="line"></div><div class="line">  <a class="code" href="classtf_1_1cudaTask.html">tf::cudaTask</a> cluster = cf.<a class="code" href="classtf_1_1cudaFlow.html#adb731be71bdd436dfb5e36e6213a9a17">kernel</a>(</div><div class="line">    (N+1024-1) / 1024, 1024, 0,</div><div class="line">    assign_clusters, d_px, d_py, N, d_mx, d_my, d_sx, d_sy, K, d_c</div><div class="line">  ).<a class="code" href="classtf_1_1cudaTask.html#ab81b4f71a44af8d61758524f0c274962">name</a>(<span class="stringliteral">&quot;cluster&quot;</span>);</div><div class="line"></div><div class="line">  <a class="code" href="classtf_1_1cudaTask.html">tf::cudaTask</a> new_centroid = cf.<a class="code" href="classtf_1_1cudaFlow.html#adb731be71bdd436dfb5e36e6213a9a17">kernel</a>(</div><div class="line">    1, K, 0,</div><div class="line">    compute_new_means, d_mx, d_my, d_sx, d_sy, d_c</div><div class="line">  ).<a class="code" href="classtf_1_1cudaTask.html#ab81b4f71a44af8d61758524f0c274962">name</a>(<span class="stringliteral">&quot;new_centroid&quot;</span>);</div><div class="line"></div><div class="line">  cluster.<a class="code" href="classtf_1_1cudaTask.html#abdd68287ec4dff4216af34d1db44d1b4">precede</a>(new_centroid)</div><div class="line">         .<a class="code" href="classtf_1_1cudaTask.html#a4a9ca1a34bac47e4c9b04eb4fb2f7775">succeed</a>(zero_c, zero_sx, zero_sy);</div><div class="line">  </div><div class="line">  <span class="comment">// we ask the executor to launch the cudaFlow by M times and then join</span></div><div class="line">  cf.<a class="code" href="classtf_1_1cudaFlow.html#a9b28ad99e4d3c0208422a2db094df277">join_n</a>(M);</div><div class="line">}).name(<span class="stringliteral">&quot;update_means&quot;</span>);</div><div class="line"></div><div class="line"><span class="comment">// ...</span></div><div class="line"></div><div class="line"><span class="comment">// build up the dependency</span></div><div class="line">h2d.<a class="code" href="classtf_1_1Task.html#a331b1b726555072e7c7d10941257f664">succeed</a>(allocate_px, allocate_py, allocate_mx, allocate_my);</div><div class="line"></div><div class="line">kmeans.<a class="code" href="classtf_1_1Task.html#a331b1b726555072e7c7d10941257f664">succeed</a>(allocate_sx, allocate_sy, allocate_c, h2d)</div><div class="line">      .<a class="code" href="classtf_1_1Task.html#a8c78c453295a553c1c016e4062da8588">precede</a>(stop);</div><div class="line"></div><div class="line">stop.<a class="code" href="classtf_1_1Task.html#a8c78c453295a553c1c016e4062da8588">precede</a>(free);</div></div><!-- fragment --><p>At the last line of the cudaFlow closure, we call <code>cf.join_n(M)</code> to ask the executor to repeatedly run the cudaFlow by <code>M</code> times and then join. Compared with the version using conditional tasking, the cudaFlow here is created only one time and thus the overhead is reduced.</p>
<div class="image">
<object type="image/svg+xml" data="kmeans_5.svg" width="100%">kmeans_5.svg</object>
</div>
<p>We can see from the above taskflow the condition task is removed.</p>
<h1><a class="anchor" id="KMeansBenchmarking"></a>
Benchmarking</h1>
<p>We run three versions of k-means, sequential CPU, parallel CPUs, and one GPU, on a machine of 6 Intel i7-8700 CPUs at 3.20GHz and a Nvidia RTX 2080 GPU using various numbers of 2D point counts and iterations.</p>
<div align="center"> <table class="markdownTable">
<tr class="markdownTableHead">
<th class="markdownTableHeadCenter">N  </th><th class="markdownTableHeadCenter">K  </th><th class="markdownTableHeadCenter">M  </th><th class="markdownTableHeadCenter">CPU Sequential  </th><th class="markdownTableHeadCenter">CPU Parallel  </th><th class="markdownTableHeadCenter">GPU (conditional taksing)  </th><th class="markdownTableHeadCenter">GPU (with predicate)   </th></tr>
<tr class="markdownTableBody" class="markdownTableRowOdd">
<td class="markdownTableBodyCenter">10  </td><td class="markdownTableBodyCenter">5  </td><td class="markdownTableBodyCenter">10  </td><td class="markdownTableBodyCenter">0.14 ms  </td><td class="markdownTableBodyCenter">77 ms  </td><td class="markdownTableBodyCenter">1 ms  </td><td class="markdownTableBodyCenter">1 ms   </td></tr>
<tr class="markdownTableBody" class="markdownTableRowEven">
<td class="markdownTableBodyCenter">100  </td><td class="markdownTableBodyCenter">10  </td><td class="markdownTableBodyCenter">100  </td><td class="markdownTableBodyCenter">0.56 ms  </td><td class="markdownTableBodyCenter">86 ms  </td><td class="markdownTableBodyCenter">7 ms  </td><td class="markdownTableBodyCenter">1 ms   </td></tr>
<tr class="markdownTableBody" class="markdownTableRowOdd">
<td class="markdownTableBodyCenter">1000  </td><td class="markdownTableBodyCenter">10  </td><td class="markdownTableBodyCenter">1000  </td><td class="markdownTableBodyCenter">10 ms  </td><td class="markdownTableBodyCenter">98 ms  </td><td class="markdownTableBodyCenter">55 ms  </td><td class="markdownTableBodyCenter">13 ms   </td></tr>
<tr class="markdownTableBody" class="markdownTableRowEven">
<td class="markdownTableBodyCenter">10000  </td><td class="markdownTableBodyCenter">10  </td><td class="markdownTableBodyCenter">10000  </td><td class="markdownTableBodyCenter">1006 ms  </td><td class="markdownTableBodyCenter">713 ms  </td><td class="markdownTableBodyCenter">458 ms  </td><td class="markdownTableBodyCenter">183 ms   </td></tr>
<tr class="markdownTableBody" class="markdownTableRowOdd">
<td class="markdownTableBodyCenter">100000  </td><td class="markdownTableBodyCenter">10  </td><td class="markdownTableBodyCenter">100000  </td><td class="markdownTableBodyCenter">102483 ms  </td><td class="markdownTableBodyCenter">49966 ms  </td><td class="markdownTableBodyCenter">7952 ms  </td><td class="markdownTableBodyCenter">4725 ms   </td></tr>
</table>
</div><p>When the number of points is larger than 10K, both parallel CPU and GPU implementations start to pick up the speed over than the sequential version. We can see using the built-in predicate of cudaFlow is two times faster than conditional tasking. </p>
</div></div><!-- contents -->
</div><!-- doc-content -->
<!-- start footer part -->
<div id="nav-path" class="navpath"><!-- id is needed for treeview function! -->
  <ul>
    <li class="navelem"><a class="el" href="Examples.html">Learning from Examples</a></li>
    <li class="footer">Generated by
    <a href="http://www.doxygen.org/index.html">
    <img class="footer" src="doxygen.png" alt="doxygen"/></a> 1.8.14 </li>
  </ul>
</div>
</body>
</html>
